Demonstration
Jean-Romain - in partnership with the Quebec ministry - has been developing approaches to detect, measure, and evaluate road detection from lidar data. His approach leverages the lidR package, which has seen massive uptake and use globally. lidR, also developed by Jean-Romain, is free and open source and is considered a key resource for lidar data processing globally.
Algorithm overview
The implementation of this algorithm is a work in progress. Function names are likely to change, methods are likely to evolve, and debug options are likely to be removed. Simply described, the current road detection and measurement algorithm measure_road() leverages the lidR processing engine to extract relevant segments of lidar data from a point-cloud with a paired linear road network. The approach is in ongoing development, but follows the following steps:
A lidar dataset and road segment and provided
The road segment is buffed and used to clip relevant lidar data
The road is divided into incremental segments
Structural lidar metrics are calculated at each segment used to define the location with the highest probability for the road
Important road measurements are calculated including drivable width and embankment width.
If the algorithm is confident that the road exists, a new (corrected) road segment is produced.
Example road
This example uses a section of SPL data and road network shapefile from the Nipissing Forest.
#--- load in required libraries ---#
library(lidR)
library(sf)
library(MFFProads)
library(dplyr)
library(leaflet)
#--- read in laz and road shapefile ---#
infolder <- "I:/Example_data/LAZ_files"
ctg <- readLAScatalog(infolder, filter = "-drop_withheld -keep_random_fraction 0.25")
roads <- st_read("E:/SPL_Roads/data/Nipissing/tst_roads_clip3.shp")
roads = st_cast(roads, "LINESTRING")mapview::mapview(list(ctg,roads), layer.name = c("Catalog","Roads"),
color = "red", map.type = "Esri.WorldImagery")Detect road and calculate metrics
Individual road segments can be visualized alongside co-located SPL data.
#--- test on single road ---#
road <- roads %>%
dplyr::filter(FID %in% c(10))Detecting the road
We define that we want to find the road using options(). An example of how data are visualized is provided below. Note the different metrics used. A spike in the purple lines denotes an increased probability of road for a given point in the transect.
#--- find the roads ---#
param = mffproads_default_parameters
param$extraction$road_buff = 25
corrected_road <- MFFProads::measure_road(road = road, ctg = ctg, dtm = dtm, water = water, param = param)mapview::mapview(list(road, corrected_road),
layer.name = c("Inaccurate", "Corrected"),
color = c("red", "blue"), map.type = "Esri.WorldImagery")
Example of road detection
Measuring the road
Following detection of the road the algorithm can also measure important road parameters.
Example of road measurement
Derive road metrics
Metrics about the road can also be measured, giving details about drivable quality and embankments. A final STATE of the roads is provided with the following gradient:
- Road exists, very confident
- Road exists, some confusion
- Road doesnt exist, some confusion
- Road doesnt exist
Example using multiple roads
This example uses a section of SPL data and road network shapefile from the Nipissing Forest.
Road segments can be visualized alongside co-located SPL data.
roads <- st_read("E:/SPL_Roads/data/Nipissing/tst_roads_clip3.shp")
roads = st_cast(roads, "LINESTRING")#--- test on single road ---#
roads <- roads %>%
dplyr::filter(FID %in% c(1:120))Detecting the roads
We define that we want to find the roads using options(). An example of how data are visualized is provided below.
#--- find the roads ---#
param = mffproads_default_parameters
param$extraction$road_buff = 25
corrected_roads <- MFFProads::measure_roads(road = roads, ctg = ctg, dtm = dtm, water = water, param = param)mapview::mapview(list(roads, corrected_roads),
layer.name = c("Inaccurate", "Corrected"),
color = c("red", "blue"), map.type = "Esri.WorldImagery")